clear
set more off
capture log close
set matsize 10000
set scheme sj

*Merge municipality data with kulturkampf data
do "prepare kulturkampf"

*Prepare interactions for acreg command

gen cath_prus = catholic_share*prussia
gen cath_oppress = catholic_share*intensity_court

*********
*Figure 3: Predicted AfD vote share in municipalities within and outside of former Prussia
*********
preserve 
keep catholic_share voters east
save "histogramdata.dta", replace
restore

preserve
acreg afd_share catholic_share prussia cath_prus i.state [weight=voters] if east == 0, latitude(_CX) longitude(_CY) dist(50) spatial


quietly margins, at(catholic_share=(0(5)100) cath_prus=(0(5)100) prussia=(0 1)) saving(predictions, replace)

use predictions, clear
list _at1 _at3 _margin in 1/5
rename _margin prediction
rename _at1 catholic_share
rename _at2 prussia
rename _at3 cath_prus

keep if cath_prus == 0 | cath_prus == catholic_share

append using "histogramdata.dta"

graph twoway (line _ci_lb catholic_share if cath_prus == 0 & prussia == 0, lcolor(gs12) lpattern(dash)) ///
(line _ci_ub catholic_share if cath_prus == 0 & prussia == 0, lcolor(gs12) lpattern(dash)) ///
(line  prediction catholic_share if cath_prus == 0 & prussia == 0, lpattern(dash)) ///
(line _ci_lb catholic_share if catholic_share == cath_prus & prussia == 1, lcolor(gs12) lstyle(connected)) ///
(line _ci_ub catholic_share if catholic_share == cath_prus & prussia == 1, lcolor(gs12) lstyle(connected)) ///
(line  prediction catholic_share if catholic_share == cath_prus & prussia == 1, lstyle(connected)) ///
(histogram catholic_share if east == 0 [fw=voters], percent bin(36) color(gs14) below), ///
graphregion(color(white)) ///
bgcolor(white) ///
yscale(range(0 15)) ///
ylabel (0(5)15) ///
xtitle("Share of Catholics") ///
ytitle("Predicted AfD Vote") ///
title("") /// 
xscale(range(0 100)) ///
legend(order(3 "Non-Prussia" 6 "Prussia"))

erase histogramdata.dta 
erase predictions.dta
restore

*********
*Figure 4: Intensity of Kulturkampf and marginal effect of Catholicism on AfD vote share
*********
keep if east == 0
reg afd_share c.catholic_share##c.intensity_court i.state [weight=voters] if east == 0, vce(cluster court_district)
margins, dydx(catholic_share) at(intensity_court=(0(1)13)) vsquish
marginsplot, xdimension(at(intensity_court)) ///
graphregion(color(white)) /// 
bgcolor(white) ///
recast(line) ///
recastci(rline) ///
ciopts(lpattern(dash))  ///
xtitle("Intensity of oppression, Court District Level") ///
ytitle("Marginal effect of Catholic") ///
yline (0) ///
title("") /// 
xscale(range(0 14)) ///
xlabel(0(2)14) ///
yscale(range(-0.2 0.2)) ///
legend(off) ///
addplot(histogram intensity_court [fw=voters] if east == 0, bin(40) color(gs14) xscale(range(0 14)) xlabel(0(2)14) fraction below)
eststo clear

*********
*Figure 6: Predicted AfD support at low and high levels of Catholicism
*********
eststo clear
reg afd_share c.catholic_share##i.prussia##c.border_rel [weight=voters] if bundesland == "Rheinland-Pfalz" & bav_prus_dist < 25 & (prussia == 1 | bavaria == 1), vce(cluster gemeindeverband)
quietly margins, at(prussia=(0 1) border_rel=(-25(0.2)25) catholic_share==(17 31 45)) vsquish saving(rdd_predictions, replace)

preserve
use rdd_predictions, replace
rename _margin prediction
rename _at1 katholisch
rename _at2 Prussia
rename _at3 border_rel

*Figure at 17% Catholic
twoway (rcap _ci_lb _ci_ub border_rel if Prussia == 1 & border_rel >0 & katholisch == 17, lcolor(gs14) recast(ci)) ///
(line  prediction border_rel if Prussia == 1 & border_rel >0 & katholisch == 17, lstyle(connected)) ///
(rcap _ci_lb _ci_ub border_rel if Prussia == 0 & border_rel <0 & katholisch == 17, lcolor(gs14) recast(ci)) ///
(line  prediction border_rel if Prussia == 0 & border_rel <0 & katholisch == 17) ///
, graphregion(color(white)) ///
bgcolor(white) ///
xline(0) ///
ytitle("Predicted AfD share") ///
yscale(range(8 18)) ///
ylabel (8(2)18) ///
xtitle("Distance to border") ///
legend(order(4 "Non-Prussia" 2 "Prussia")) ///
title("Predicted at 17 percent Catholic population", ring(61) size(medlarge)) saving(fig6a, replace)


*Figure at 45% Catholic
twoway (rcap _ci_lb _ci_ub border_rel if Prussia == 1 & border_rel >0 & katholisch == 45, lcolor(gs14) recast(ci)) ///
 (line  prediction border_rel if Prussia == 1 & border_rel >0 & katholisch == 45, lstyle(connected)) ///
(rcap _ci_lb _ci_ub border_rel if Prussia == 0 & border_rel <0 & katholisch == 45, lcolor(gs14) recast(ci)) ///
(line  prediction border_rel if Prussia == 0 & border_rel <0 & katholisch == 45) ///
, graphregion(color(white)) ///
bgcolor(white) ///
xline(0) ///
ytitle("Predicted AfD share") ///
yscale(range(8 18)) ///
ylabel (8(2)18) ///
xtitle("Distance to border") ///
legend(order(4 "Non-Prussia" 2 "Prussia")) ///
title("Predicted at 45 percent Catholic population", ring(61) size(medlarge)) saving(fig6b, replace)


graph combine fig6a.gph fig6b.gph, commonscheme ///
	graphregion(c(white) lc(white)) plotregion(m(medlarge)) scheme(s2mono) ysize(6) xsize(14)
restore

erase fig6a.gph
erase fig6b.gph

*********
*Figure 7: Regional strength of the Volksverein and marginal effect of Catholicism on AfD vote
*********
*Figure Volksverein 1913
reg afd_share c.catholic_share##c.volksverein_1913 i.state [weight=voters] if east == 0, vce(cluster diocese)
margins, dydx(catholic_share) at(volksverein_1913=(1/7)) vsquish
marginsplot, xdimension(at(volksverein_1913)) ///
graphregion(color(white)) /// 
bgcolor(white) ///
recast(line) ///
recastci(rline) ///
ciopts(lpattern(dash))  ///
xtitle("Strength of Volksverein 1913, Diocese Level") ///
ytitle("Marginal effect of Catholicism") ///
yline (0) ///
title("") /// 
xscale(range(0 7)) ///
xlabel(0(1)7) ///
yscale(range(-0.1 0.25)) ///
ylabel(-0.1(0.1)0.2) ///
legend(off) ///
addplot(histogram volksverein_1913 [fw=voters], bin(40) color(gs14) xscale(range(0 7)) xlabel(0(1)7) yscale(range(-0.13 0.25)) ylabel(-0.1(0.1)0.2) fraction below) saving(fig7a, replace)

*Figure Volksverein 1927
reg afd_share c.catholic_share##c.volksverein_1927 i.state [weight=voters] if east == 0, vce(cluster diocese)
margins, dydx(catholic_share) at(volksverein_1927=(0(0.2)4.2)) vsquish
marginsplot, xdimension(at(volksverein_1927)) ///
graphregion(color(white)) /// 
bgcolor(white) ///
recast(line) ///
recastci(rline) ///
ciopts(lpattern(dash))  ///
xtitle("Strength of Volksverein 1927, Diocese Level") ///
ytitle("Marginal effect of Catholicism") ///
yline (0) ///
title("") /// 
xscale(range(0 4.5)) ///
xlabel(0(1)4) ///
yscale(range(-0.1 0.25)) ///
ylabel(-0.1(0.1)0.2) ///
legend(off) ///
addplot(histogram volksverein_1927 [fw=voters], bin(40) color(gs14) xscale(range(0 4.5)) yscale(range(-0.13 0.25)) ylabel(-0.1(0.1)0.2) fraction below) saving(fig7b, replace)


graph combine fig7a.gph fig7b.gph, commonscheme ///
	graphregion(c(white) lc(white)) plotregion(m(medlarge)) scheme(s2mono) ysize(6) xsize(14)

erase fig7a.gph
erase fig7b.gph	
	
*********
*Figure 8: Participation at Katholikentag festivals and marginal effect of Catholicism on AfD share
*********
reg afd_share c.catholic_share##c.katholikentag_est i.state [weight=voters] if east == 0, vce(cluster diocese_1990)
margins, dydx(catholic_share) at(katholikentag_est=(-1(0.1)1.5)) vsquish
marginsplot, xdimension(at(katholikentag_est)) ///
graphregion(color(white)) /// 
bgcolor(white) ///
recast(line) ///
recastci(rline) ///
ciopts(lpattern(dash))  ///
xtitle("Participation at Katholikentage") ///
ytitle("Marginal effect of Catholicism") ///
yline (0) ///
title("") /// 
xscale(range(-1 1.5)) ///
xlabel (-1(0.5)1.5) ///
legend(off) ///
addplot(histogram katholikentag_est [fw=voters] if east == 0, bin(25) color(gs14) xscale(range(-1 1.5)) fraction below)

*********
*Figure 9: Centre Party vote shares in Reichstag elections, West Germany
*********
clear
use "dataset_german_empire.dta", replace

foreach x in 1871 1874 1877 1878 1881 1884 1887 1890 1893 1898 1903 1907 1912 {
          quietly eststo r`x': reg `y' zentrumshare c.catholic_share##i.prussia agr_empl_1895 ind_empl_1895 urban_pop_1900 landgini1895 non_german_sp_1900  if yr == `x'  & east == 0
        }	
		
coefplot (r187*, asequation \ r188*, asequation \ r189*, asequation \ r19*, asequation) , drop(_cons agr_empl_1895 ind_empl_1895 urban_pop_1900 landgini1895 non_german_sp_1900) xline(0) swapnames ///
graphregion(color(white)) bgcolor(white) ///
coeflabels(r1871 = "1871" r1874 = "1874" r1877 = "1877" r1878 = "1878" r1881 = "1881" r1884 = "1884" r1887 = "1887" r1890 = "1890" ///
r1893 = "1893" r1898 = "1898" r1903 = "1903" r1907 = "1907" r1912 = "1912", labsize(vsmall)) ///
rename(catholic_share = "Percent Catholic" 1.prussia = Prussia 1.prussia#c.catholic_share = "Prussia x Catholic") ///
title("Regression of centre party vote share on Prussia and Catholicism, 1871-1912", size(medsmall) margin(b=5))


*********
*Figure 10: Share of West Germans most opposed to the AfD
*********
clear
use "dataset_individual.dta", replace

gen AfD_resist = thermometer_AfD_bin - 1
replace AfD_resist = (-1)*AfD_resist

preserve
keep if eastwest == 1
cibar AfD_resist, over(prussia_bin catholic) level(83) graphopts(yscale(range (0 0.8))ylabel(0(0.2)0.8) ///
ytitle(Share of most opposed to AfD) graphregion(color(white)) bgcolor(white)) ///
barlabel(on) blf(%9.3f) blgap(0.05)
restore
